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and surface volcanic deposits, which has remained unexplained since its discovery in the Apollo era. The impact- 
induced heat drives hemisphere-scale mantle convection, which would sequester Th- and Ti-rich lunar magma 
ocean cumulates in the nearside hemisphere within a few hundred million years if they remain immediately 
beneath the lunar crust at the time of the SPA impact. A warm initial upper mantle facilitates generation of a 
pronounced compositional asymmetry consistent with the observed lunar asymmetry. 


INTRODUCTION 

One of the fundamental mysteries of lunar evolution has been the 
origin of the prominent and ancient asymmetry between the Moon’s 
farside and nearside hemispheres (Fig. 1) (1). The giant South 
Pole-Aitken (SPA) lunar impact basin overprints the farside 
hemisphere’s generally high albedo—diagnostic of the anorthositic 
lunar primary crust—and low abundances of elements such as Th 
and Ti. In contrast, Th and Ti abundances are elevated across the 
nearside, and swaths of low albedo covering a similar area denote 
lunar maria—dense basaltic lavas that preferentially resurfaced the 
nearside beginning ~4.3 billion years ago (Gya) (2-4). The nearside’s 
most extensive maria and highest Th concentrations encompass the 
SPA basin antipode and extend toward the southwest, marking the 
Procellarum KREEP Terrane (PKT) (2), a province that hosts a 
unique geochemical component with elevated levels of potassium 
(K), rare earth elements (REE), phosphorous (P), and other trace 
and heat-producing elements such as Th (5, 6). 

KREEP originated from the final magma ocean-derived residu- 
um that crystallized beneath the lunar crust, alongside other late- 
stage cumulates bearing Ti-oxide minerals like ilmenite (5-7). The 
lunar compositional asymmetry is commonly hypothesized (1, 8-11) 
to have originated through the redistribution of a global layer of 
subcrustal late-stage cumulates that occurred as early as ~4.4 Gya 
(the end of lunar magma ocean crystallization (5-7)). However, the 
catalyst for such a redistribution has remained enigmatic. Here, we 
show that the impact event that formed the SPA basin by 4.3 Gya 
(12) was appropriately timed and imparted sufficient heating to 
initiate lunar mantle dynamics that entrain subcrustal KREEP- and 
Ti-rich late-stage cumulates in lateral flow toward the nearside 
hemisphere. Under warm upper mantle conditions favored by SPA 
basin formation constraints (13), a pronounced compositional 
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asymmetry forms that is consistent with surface observations of the 
ancient lunar asymmetry. 

We examine lunar evolution under two end-member scenarios 
for the state of the interior at the time of the SPA impact (Fig. 2, 
first row). All numerical convection simulations in this study 
(see Supplementary Text) begin with a global layer of ilmenite-bearing 
late-stage cumulates overlying the lunar mantle at the time of the 
SPA impact, in accordance with recent remote sensing analyses 
(14). Reference scenario 1 (RS1) (Fig. 2A) considers an SPA impact 
into the stratified mantle predicted to exist immediately following 
solidification of the primordial lunar magma ocean (7). On the 
basis of experimental constraints on the rheology of ilmenite-bearing 
cumulates (15), the unmixed mantle of RS1 allows for late-stage 
cumulates that are weakened (i.e., lower viscosity) by a moderate 
ilmenite content (see Materials and Methods). Alternatively, for 
RS2 (Fig. 2, B and C), we consider an SPA impact into a thermally 
homogenized post-magma ocean cumulate mantle, in line with 
constraints for the lunar mantle thermal structure at the time of 
impact [e.g., (13)]. To account for the possibility of some ilmenite 
depletion by mixing of the underlying mantle (15), the subcrustal 
late-stage cumulate layer is less ilmenite-rich than in RS1 and has 
no compositional weakening effect. However, despite the difference 
in compositional effects, the temperature dependence of viscosity in 
our simulations means that the warm upper mantle of RS2, including 
the late-stage cumulate layer, flows more easily than the cooler 
upper mantle of RS1. 

RS1 and RS2 consider best-fit SPA impact parameters from 
vertical basin formation simulations (85-km radius projectile, 10-km/s 
impact) (13, 16). The more likely formation of the elliptical SPA 
basin by an oblique impact (17) allows for a more energetic collision 
[e.g., (18)] with more extensive heating, because basin radius is a 
function of the vertical component of impact velocity (19). To 
further highlight the range of possible outcomes, RS2B allows for 
an approximately order of magnitude more energetic SPA impact 
with the same initial mantle conditions as RS2 (Fig. 2). 


RESULTS AND DISCUSSION 
At the start of our simulations, we determine impact shock-induced 
heating using an analytic model that scales with the energy of the 
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Fig. 1. Maps showing the lunar surface compositional asymmetry. (A) Lunar 
Reconnaissance Orbiter Wide Angle Camera mosaic (39). (B) Surface Th concentration 
(40) [used as a proxy for KREEP distribution (2)]. (C) Surface Ti concentration (47). 
Left: Lunar nearside hemisphere. Right: Farside hemisphere. South Pole—Aitken 
(SPA) basin ellipse (dashed white line) and basin center antipode (cyan circle) (77) 
are marked. Procellarum KREEP Terrane (PKT) outline (pink) is denoted by the 
3.5-parts per million Th contour (2). Maps are plotted in Lambert equal-area 
azimuthal projection. 


impact (20, 21). We find, for the range of plausible mantle pa- 
rameters, that the SPA impact emplaces a mantle thermal anomaly 
(Fig. 2, second row) that exerts control over lunar interior evolution 
for hundreds of millions of years. The SPA-induced thermal anomaly 
drives a hemisphere-scale upwelling beneath the SPA basin on the 
lunar farside (Fig. 2, third and bottom rows), and the larger thermal 
anomaly of RS2B leads to a similarly sized downwelling beneath the 
impact antipode within the nearside PKT region (location marked 
in Fig. 1). 

In every scenario examined, the SPA-induced hemisphere-scale 
convection pattern brings about lateral upper mantle flow toward 
the SPA antipode that entrains late-stage cumulates and generates a 
hemispheric compositional asymmetry, which is especially pronounced 
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given warmer initial upper mantle conditions (Fig. 3). The de- 
velopment of the asymmetry in KREEP- and ilmenite-bearing 
late-stage cumulates as a result of the SPA impact agrees with the 
timing and distribution of multiple observed features of the lunar 
asymmetry. 

An early, global, subcrustal late-stage cumulate layer is supported 
by evidence for a global lunar magma ocean and the confirmed 
presence of Ti-rich basalts and KREEP on the lunar nearside (5-7) 
as well as recent remote sensing analyses, which indicate KREEPy 
late-stage cumulate lithologies that are associated with Th anomalies 
in the SPA basin (14). The expected lateral and vertical distributions 
of emplaced SPA ejecta [e.g., (22)] are proposed by Moriarty et al. 
(14) to agree with the distributions of observed lithologies and surface 
Th on the farside after processes including postimpact basin modi- 
fication and inward slumping of Th-rich ejecta, thus supporting the 
presence of late-stage cumulates that were intact beneath the lunar 
crust and excavated during formation of the SPA basin. 

In addition, crater retention ages suggest that the PKT is younger 
than the SPA basin (12), and the most ancient nearside mare basalts 
erupted as early as 4.3 Gya (3, 4) [up to ~200 million years (Ma) 
after the proposed age of the SPA impact (12)] followed by intense 
episodes of nearside mare volcanism beginning ~3.8 Gya (23) 
[~500 Ma to 700 Ma after the SPA impact (12)]. The timing and 
location of the intense nearside mare volcanism are well explained 
by an early sequestration of heat-producing KREEP in the lunar 
nearside (1). Our findings show that within 300 Ma (Fig. 2, third 
row) and until at least 600 Ma after the SPA impact (Figs. 2, bottom 
row, and 3), late-stage cumulates are concentrated in the nearside, 
which would lead to a KREEP and thermal asymmetry that could be 
linked with the formation of the PKT and nearside mare basalts. 

Toward the end of our simulations, at ~3.9 Gya (400 Ma to 
600 Ma after the SPA impact) (12), the Imbrium impact would have 
excavated the locally concentrated KREEP material on the nearside 
(Fig. 3). Localization of KREEP-rich material to the lunar nearside 
followed by its excavation during the Imbrium impact has previously 
been suggested to have generated the Th signature across much of 
the lunar surface (24). 

The simulations shown in Fig. 2 represent outcomes along a 
continuum where the degree of nearside KREEP localization 
catalyzed by the SPA basin-forming impact is dependent on the 
thermal state and rheology of the lunar mantle as well as the energy 
of the SPA impact. For the relatively cooler and more rigid near- 
surface conditions of the stratified mantle in RS1, lateral transport 
of late-stage cumulates is moderate (farside thickness is reduced to 
about half and nearside thickness increases slightly), and the layer 
remains intact beneath the crust (Figs. 2 and 3). However, for 
warmer and less rigid upper mantle conditions predicted at the time 
of the SPA impact [e.g., (13)], as in RS2 and RS2B, nearly all late- 
stage cumulates are transported away from the farside and accumulate 
concentrically around the SPA antipode. Constraints on the origin 
of Th-rich ejecta on the lunar surface suggest that basin-forming 
impacts other than Imbrium (i.e., those outside of the PKT) exca- 
vated little to no KREEP (24), and the majority of extensive maria 
are limited to the area within and proximal to the PKT (Fig. 1); 
these observations point to the pronounced late-stage cumulate 
localization demonstrated by RS2 and RS2B. 

With a warm upper mantle as well as a more energetic SPA 
impact than constrained by axisymmetric models (13, 16), late-stage 
cumulates may have been later driven deep into the lunar mantle. 
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Fig. 2. Cross sections of 3D lunar interior evolution simulations showing the effects of the SPA impact. All cross sections are oriented with the farside SPA basin 
at the north pole. The lunar core (inner white circle) is shown to scale, while the crust (dark gray outer layer) is illustrated in postprocessing with exaggerated thickness. 
(A) Reference scenario 1 (RS1) (stratified post-magma ocean mantle) shown immediately before (top row) and after (second row) the SPA impact as well as 300 Ma (third row) 
and 600 Ma (bottom row) after impact. Left: Composition, distinguishing for simplicity between the material that is primarily late-stage cumulates (pink) and primarily 
other lunar mantle (gray). The dashed circle denotes ~500-km depth, the approximate maximum source depth of Ti-rich igneous rocks taken from the lunar nearside 
surface (42). Right: Temperature. (B) Same as (A) for RS2 (mixed mantle after magma ocean solidification). (C) Same as (A) for RS2B (mixed mantle, more energetic impact). 
(D) Schematic cross-sectional illustrations of SPA-induced shock heating and simplified SPA-induced convection pattern. Shock heating at the time of the SPA impact 
(second row) is represented with an arbitrary isotherm (dashed orange line), radially symmetric from a point at depth below the impact. 


An SPA-catalyzed downwelling of late-stage cumulates as highlighted 
by RS2B (Fig. 2C) would be consistent with several lines of evidence 
that suggest that dense, subcrustal, KREEP- and ilmenite-bearing 
lithologies were transferred deeper into the lunar interior. Volcanic 
glasses found on the lunar nearside—with melt source depths con- 
strained to ~500-km—assimilated high-Ti and KREEP compo- 
nents, necessitating the sinking of late-stage cumulates (25). 
Reanalysis of lunar seismic data indicates a low-viscosity layer sur- 
rounding the lunar core (26), which is proposed to be consistent 
with a stable layer of ilmenite-bearing lunar magma ocean cumu- 
lates that foundered to the core-mantle boundary (CMB) (15, 27). 
Furthermore, constraints from crater retention ages and remote sensing 
analyses suggest that the SPA impact occurred before the removal of 
early late-stage cumulates from beneath the crust (12, 14), allowing 
for an SPA-induced downwelling of late-stage cumulates. 

A nearside downwelling of dense ilmenite-bearing cumulates, 
due to passive gravitational instability, has previously been put 
forth to explain the origin of the lunar asymmetry (8), although the 
development of such a hemisphere-scale downwelling depends 
greatly on the rheology of late-stage cumulates (8, 15). A passive 
density-driven mechanism for the downwelling of late-stage cumu- 
lates has also been demonstrated as unlikely to entrain the KREEP-rich 
magma ocean residuum (28). Our results show that hemisphere- 
scale convection in the aftermath of the SPA impact (a known 
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basin-forming event) can initiate downwelling of late-stage cumulates 
even if the layer is not weaker or denser than the mantle. SPA- 
induced convection therefore offers a viable stand-alone catalyst 
that could be augmented by passive gravitational instability. 

Our simulations of interior evolution for plausible early lunar 
mantle states demonstrate that heating from the SPA basin-forming 
impact initiates long-lived hemisphere-scale mantle convection, 
which, especially under warm mantle conditions favored by SPA 
basin formation constraints (13), is capable of generating a compo- 
sitional asymmetry consistent with noted observations. Although 
the major dynamical influence of interior heating by large impacts 
has already been examined extensively for Mars [e.g., (29, 30)], 
near-surface processes discussed here may provide insight to the 
giant impact origin hypothesis of Mars’s north-south hemispheric di- 
chotomy [e.g., (31, 32)]. Furthermore, neither the reorganization of 
lunar interior dynamics nor the nearside concentration of late-stage 
cumulates caused by SPA impact heating have been widely con- 
sidered in investigations of lunar evolution. Our results presented 
here indicate that SPA-induced convection is a fundamental con- 
sideration for lunar history, and the associated hemisphere-scale 
processes could have many implications beyond what is discussed here. 
For example, the nearside localization and possible downwelling of 
dense late-stage cumulates may influence the offset between the 
Moon’s center of figure and center of mass. 
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Fig. 3. Plots of the simulated and observed lunar compositional asymmetry. (Left axis) Simulation results for surface-equivalent thickness of KREEP- and Ti-bearing 
lunar late-stage cumulates (pink) at 600 Ma after the SPA impact, azimuthally averaged with respect to distance from the SPA basin. (Right axis) Mean (black line) and 16 
variation (gray shaded region) of present-day lunar surface Th concentration (40) azimuthally averaged with respect to distance from the SPA basin. Simulated late-stage 
cumulate thicknesses are normalized to the initial layer thickness and are calculated from the total amount of late-stage cumulates above 500-km depth, the approximate 
maximum source depth of nearside Ti-rich igneous rocks (42). The left and right vertical axes are independently scaled to show the full range of each dataset. Extent of 
the SPA basin (77) and approximate extent of the Th-rich PKT are denoted as blue shaded regions. The Imbrium (IMB) basin center position (43) is marked by a vertical 


dashed line at the shortest great circle distance from the center of SPA. 


MATERIALS AND METHODS 

Thermochemical evolution 

Our lunar interior evolution simulations use a modified version of 
CitcomS, a three-dimensional, finite-element, thermochemical evo- 
lution code with Lagrangian tracer—based composition tracking 
(33-35). We consider a lunar mantle that can be described using the 
Boussinesq approximation and infinite Prandtl number. The govern- 
ing equations of conservation of mass, momentum, energy, and 
composition take the following nondimensional forms (33, 34) 


uj, =0 
-P; + (n(ui, + uj); + RaT8ir =0 
T,+ujTj= Ti; + Q 


Cr+ uj;C;= 0 


where u is velocity, P is dynamic pressure, n is viscosity, Ra is the 
Rayleigh number, T is temperature, 5;, is the Kronecker delta (where 
subscript r denotes the radial direction), t is time, Q is the internal 
heat production, and Cis composition. Subscripts i and j are spatial 
indices, and a subscript comma indicates the partial derivative with 
respect to the following index or variable. Reference values used in 
nondimensionalization are listed in table S1. 

In CitcomS, the basal heating Rayleigh number Ra is de- 
fined as (33) 
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where AT is the temperature drop from the CMB to the surface (in units 
of K), Ro is the planetary radius, and subscript 0 denotes the dimen- 
sional reference value of a quantity used in nondimensionalization. 
Relative to typical definitions of Ra, which normalize with the con- 
vective layer thickness, the CitcomS definition using the full radius 
Ro leads to a value of Ra that is approximately 20% larger for the 
Moon given our selected lunar parameters (table S1). As detailed in 
the Lunar interior models section, our simulations exclude the lunar 
crust, so AT and Ro are defined so that the outer surface of the model 
is the boundary between the crust and the underlying mantle 
(including late-stage cumulates). Hence, the model domain is the 
extent of the lunar mantle, from the CMB to the crust-mantle 
boundary. Boundary conditions at both surfaces are isothermal and 
free slip. The effects of rigid body rotation are removed. 

Our simulations assume Newtonian rheology where temperature-, 
depth-, and composition-dependent viscosity is computed according 
to the dimensionless equation 


E, + Va(1 — r) 
T+T. 


Tin) = Argesn( E, + Va(1 - | 


surf 1+ surf 


where A is a composition-dependent viscosity prefactor, 14 is 
the dimensionless reference viscosity (with a magnitude of unity 
under the nondimensionalization scheme of CitcomS), E, is the 
activation energy, V, is the activation volume, T.,,. is the dimen- 
sionless temperature at the outer surface of the model, R, is the core 
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radius, and all other symbols are as defined above. Dimensionless 
activation energy E, and volume V, are defined as 


where E,, and V,, are the dimensional reference values for activa- 
tion energy and activation volume, respectively, R is the universal 
gas constant, and all other symbols are as defined above. 

The domain of our models consists of a spherical shell divided 
into 12 spherical caps, each with a grid resolution of 21 by 21 by 
39 nodes in the directions of latitude, longitude, and radius, respec- 
tively. The global mesh is composed of (20 x 20 x 38) x 12 = 182,400 
brick elements with radial length of ~36.1 km and lateral length 
ranging from ~18.3 km at the CMB to ~94.4 km at the outer surface 
of the model. 

Composition is tracked using advective tracers, which allow for 
composition-dependent viscosity (34). Our simulations use two 
tracer compositions—one for late-stage cumulates and one for the 
rest of the lunar mantle. Simulations are initialized with 50 tracers 
per element. At each time step of a simulation, the ratio of each 
tracer composition is used to compute the viscosity reduction 
coefficient A in each brick element of the mesh as a geometric mean 
between the late-stage cumulate compositional viscosity prefactor 
Ajsc and unity. 

Time step size is computed such that the size of each step is 70% 
of the dynamically computed maximum stable advection time step. 
For our model parameters, this yields time steps on the order of ~10"' 
Ma to 10° Ma. 


Impact-induced shock heating 

Large impacts like that which formed the SPA basin impart sub- 
stantial instantaneous heating to planetary interiors, and the resultant 
thermal anomaly can significantly influence subsequent mantle dy- 
namics (20, 21). We model the emplacement of the SPA impact- 
induced thermal anomaly (impact parameters listed in table S1) 
by prescribing a temperature increase throughout the lunar man- 
tle according to the Hugoniot release method of Watters et al. (21) 


aTsPs)= 2( x85 f?) (S)'v Inf 1)) 


where 


flPs)= 


-1 
2SP5 (: ASPs 1) 
C“ po C" po 
and where shock-induced temperature increase AT; is a function of 
shock-induced pressure increase Ps, cp is the specific heat capacity, 
Po is the reference density, and S and C are the slope and intercept 
of the linear Hugoniot shock equation of state (EOS), respectively 
(Cis also the acoustic velocity in the target material). 
Instantaneous shock pressure Ps (hence AT;) is approximately 
uniform within a spherical region referred to as the isobaric core, 
outside of which Ps and AT, decay exponentially with radial distance 
(20, 21). The isobaric core is located at depth d- and has radius r;, 
(20). For the impact parameters used in this study (table S1) 
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where values for the radius of the impactor Rimp = 85 km and velocity 
of the impact Vimp = 10 km/s are from best-fit constraints of SPA 
basin formation simulations (13, 16). The pressure within the 
isobaric core, P,, is 


P. = Po(C + Suc) ue = 219.5 GPa 


where u, is the shock-induced particle velocity within the isobaric 
COre Ue = +(Vimp) under the assumption that the target and projectile 
materials are the same, e.g., dunite (21), and other symbols are as 
defined above. Outside this region, the peak shock pressure, Ps, 
decays according to the inverse power law 


Pir) = Pe (*) _ PoC + Suc) Uc ($) al 


where r is the radial distance from the origin of the isobaric core, 
and other symbols are as defined above. 

P3 at a given point is defined as the difference between the peak 
shock pressure P; (or P, when r < r,) at that point and the preimpact 
pressure Pp at that point, ie., Ps = Ps — Po. 

Our modified version of CitcomS calculates shock-induced 
temperature increase at each node of the finite element mesh at an 
impact time defined relative to the simulation start time (our simu- 
lations use timpact = 5 Ma after initialization). Once the start time 
(units of Ma) of a step exceeds the specified impact time and the 
energy equation has been solved for that step, shock heating is 
calculated and added to the temperature at each node of the mesh. 
Following the example of Roberts and Arkani-Hamed (36), we 
assume that any melt produced in the impact returns to the mantle 
solidus within one time step because the time scale of crystallization 
of a body of impact-generated melt at the surface is much quicker 
than the time scale of solid-state convection (37). That is, the maximum 
temperature Tmax Calculated during our shock heating routine at 
any node is limited by the mantle solidus. 

From the equations and assumptions defined above, the magni- 
tude and extent of the lunar mantle thermal anomaly generated by 
the SPA impact depend in part on the lunar mantle shock EOS 
(from S and C) and the energy of the SPA impact (from Rimp and 
Vimp). Across a range of materials potentially relevant to the lunar 
mantle, differences in constrained values of S and C lead to little 
enough variation in the SPA-induced thermal anomaly (fig. $1) 
that we expect similar outcomes with any reasonable choice of 
shock EOS. To produce representative scenarios (Figs. 2 and 3), 
we therefore select values that lead to an intermediate thermal 
anomaly within the range of potentially relevant materials [S = 1.5 
and C= 5.2 km/s (38)]. 


Lunar interior models 

To examine the first-order effects of the SPA-induced thermal 
anomaly on lunar interior dynamics, focusing especially on the 
evolution of an initially global, subcrustal layer (5-7, 14) of KREEP- 
and Ti-bearing late-stage cumulates, we consider a two-layer lunar 
mantle (fig. $2). At the start of our models, a 60-km-thick layer of 
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late-stage cumulates overlies a compositionally uniform silicate 
mantle that extends to the lunar CMB. Because of the temperature 
dependence of viscosity in our models and the expected low tem- 
perature of the lunar crust (7, 13), the crust is taken to be much 
more viscous than the mantle and is excluded from our models. 
This also reduces the overall viscosity contrast in the simulations 
and ensures that CitcomS is able to maintain numerical stability. 
Compositional density is uniform throughout the model. A dense 
late-stage cumulate layer could alter our simulations by promoting 
Rayleigh-Taylor instability of the layer [e.g., (8, 28)]; however, we 
expect our conclusions to be largely unaltered as short-wavelength 
instabilities would likely be entrained in the stronger SPA-induced 
lateral flow field and long-wavelength instabilities may augment the 
compositional redistribution that we observe and facilitate a near- 
side downwelling of late-stage cumulates. 

Simulations are initialized with one of two end-member thermal 
structures (fig. $2) that represent the range of plausible lunar mantle 
thermal states at the time of the SPA impact. Mantle model 1 allows 
for an SPA impact that occurs immediately after solidification of 
the lunar magma ocean and formation of the late-stage cumulate 
layer beneath the crust. The thermal structure of model 1 is adapted 
from model results for the stratified lunar mantle state at the end of 
magma ocean solidification (7). Mantle model 2 considers a 
thermally homogenized cumulate portion of the post-magma ocean 
lunar mantle, in line with some considerations of the lunar thermal 
state at the time of the SPA impact [e.g., (13)]. The late-stage cumu- 
late layer is initially intact beneath the lunar crust in both models on 
the basis of recent remote sensing analyses that suggest SPA exca- 
vated a subcrustal KREEP-rich layer (14). The deep, early lunar 
thermal structure is not well constrained; our models allow for a 
warm, basal primordial mantle layer that remains unmixed with the 
overlying magma ocean-derived cumulate mantle at the time of 
the SPA impact, based on models of magma ocean solidification 
(7), but thermal constraints relevant for the upper mantle (13) 
simply assume an approximately adiabatic profile all the way to 
the CMB. 

Since the temperature of the outer surface of each mantle model 
is different, AT also differs for each (table S2). The basal heating 
Rayleigh number Ra and nondimensional heating number Q thus 
differ between each, where 


HR; 
Q=—aF 
KcpAT 


and where H is internal heating in terms of mass. For mantle model 1, 
Ra = 1.3 x 10° and Q = 13.72. For mantle model 2, Ra = 4.4 x 10° and 
Q=41.16. The influence of internal heating and the convective driving 
force of the SPA-induced thermal anomaly mean that Ra as reported 
here underrepresents the vigor of convection. 

We additionally consider compositional viscosity reductions on 
the late-stage cumulate layer. Recent rheological experiments indicate 
that ilmenite fractions of approximately 5 to 15%, as expected for 
lunar late-stage magma ocean cumulates, are capable of reducing 
the viscosity of silicate rock by up to ~2 to 3 orders of magnitude 
due to the low viscosity of ilmenite (15). To broadly characterize the 
relationship between relative ilmenite content and the outcome of 
SPA-induced convection, we examine model scenarios with zero-, 
one-, and two-order magnitude viscosity reductions on the late- 
stage cumulate layer. 
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